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In this talk, we present our recent work on next-to-leading order (NLO) nuclear parton 
distribution functions (nPDFs), which we call EPS09. As an extension to earlier NLO 
analyses, we complement the deep inelastic scattering and Drell-Yan dilepton data by 
inclusive midrapidity pion production measurements from RHIC to reduce the other- 
wize large freedom of the nuclear gluon densities. In addition, our Hessian-type error 
analysis leading to a collection of nPDF error sets, is the first of its kind among the 
nPDF analyses. 

1 Introduction 

The global analyses of the free nuclcon parton distribution functions (PDFs) are based on 
the asymptotic freedom of QCD, parton evolution and factorization. These features provide 
the justification for computing hard-processes cross-sections, schematically as 

where fiS denote the scale-dependent PDFs, and &ij^h+x are perturbatively computable 
coefficients. The factorization theorem has turned out to work extremely well with more 
and more different types of data included in the free proton analyses. In the case of bound 
nucleons factorization is not as well-established, but it has nevertheless proven to provide 
a very good description [U [51 [31 [H O [5] of the observed nuclear modifications crbound/cfree 
in deep inelastic scattering (DIS) and Drell-Yan (DY) dilepton production involving nuclear 
targets. Here, we summarize our recent NLO analysis of the nuclear PDFs and, in particular, 
their uncertainties 

2 Analysis method and framework 

Our analysis follows very similar pattern as most free proton analyses do: 

A. The PDFs are parametrized at an initial scale Qq imposing the sum rules. 

In this work we do not parametrize the absolute PDFs, but rather the x and A dependences 
of the nuclear modification factors Rf{x, Qq) on top of a fixed set of free proton PDFs: 

/^(x,Q^)^i?f(.x,Q^)/™i^(x,02). 

Above f^'^^'^^-^^{x,Q^) refers to a CTEQ set of the free proton PDFs [8] m the zero- 
mass variable flavour number scheme, and we consider three different modification factors: 
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Ryix, Qq) for both u and d valence quarks, Rg{x, Qq) for all sea quarks, and Rq{x, Qq) for 
gluons. 

B. The nuclear PDFs are evolved to other perturbative scales > Qq by 
the DGLAP equations. An cfScient numerical solver for the parton evolution is an 
indispensable ingredient for any parton analysis, but in the case of nuclear PDFs this is even 
more critical as we need to repeatedly do the evolution separately for 13 different nuclei. On 
the other hand, the relatively low Q^-values spanned by the data utilized, < 200 GeV^, 
and the fact that they come as ratios crbound/ffrcc which are more stable against small 
evolution inaccuracies, make such task somewhat easier. Our DGLAP code is based on a 
semi-analytical method described e.g. in [TU]. 

C. The cross-sections are computed using the factorization theorem. 

D. The computed cross-sections are compared with the experimental data, 
and the initial parametrization is varied to establish an optimal fit to the data. 

We define the best agreement as the minimum of a generalized x^-function 

Within each data set A, the Di denotes the experimental data value with cr^ point-to-point 
uncertainty, and Ti is the theory prediction corresponding to a parameter set {a}. For the 
pion data, the PHENIX experiment estimates an overall ^ 10% normalization uncertainty 
ct5V°™, and the normalization factor S [1 - o-]^°"°, 1 + 0-^°™] is non-trivial, i.e. 7^ 1. 
Its value is determined by minimizing Xn and the final Jn is thus an output of the analysis. 
The weight factors wpj are used to amplify the importance of those data sets whose content 
is physically relevant, but contribution to would otherwize be too small to be noticed by 
an automated minimization. 

E. The uncertainties are estimated. Besides finding the central set of parameters 
{aP} that optimally fits the data, quantifying the uncertainties stemming from the experi- 
mental errors has become an integral part of the modern PDF fits. In this work, we employ 
the Hessian method |11| , which rests on a quadratic approximation 

ij * ^ ij 

around the vicinity of the minimum Xo- Non-zero off-diagonal elements in the Hessian 
matrix H defined above are a signal of correlations between the original fit parameters and 
it is useful to make a change of variables {a} — > {z} that diagonalizcs the Hessian matrix. 
Constructing the so-called PDF error sets 5*^ is what ultimately makes the Hessian method 
so useful. Each is obtained by displacing the fit parameters to the positive/negative 
direction along the eigenvector Zk of the Hessian matrix such that x^ grows by a certain 
amount Ax^- Using these sets, the upper and lower uncertainty of a quantity X can be 
written e.g. as 

{AX+r « ^[max{A(5+)-A(5"),A(5,7)-A(5°),0}]' (1) 

k 

(AX-)' « ^[max{A(^°)-A(5+),A(5°)-A(5,7),0}]', 

k 



<Ji 
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where X{S^) 



denotes the value of the quantity X computed by the set Sf and where S'^ 



is the best fit. Requiring each data set to remain close to its 90%-confidence range, we end 
up with a choicc0 Ax^ = 50. 
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3 Results and Conclusions 



We briefly go through 
the main results from the 
present analysis, starting 
with Fig. [T] where we 
plot the obtained mod- 
ifications for Lead at 
two scales making their 
scale-dependence thereby 
clearly visible. It should 
be noticed that even a 
rather large uncertainty 
band at small- a; gluons 
shrinks in the scale evolu- 
tion — a clear prediction 
of the DGLAP approach 
that might be testable in 
the future colliders. 

As the DIS and DY 
data constitute the major 
part of the available ex- 
perimental data we dis- 
play in Fig.[2]some repre- 
sentative examples of the 
measured nuclear modi- 
fications with respect to 
Deuterium, 
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Figure 1: The obtained nuclear modifications for Lead at the 
initial scale = 1.69 GeV^ and at lOOGeV^. The thick 

black lines indicate the best-fit, whereas the dotted green curves 
denote the individual error sets which combine to the shaded 
bands like in Ec^. ([1]). 
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for different nuclei compared with the EPS09. The shaded blue bands denote the uncer- 
tainty propagated from the 30 EPS09 error sets and, as should be emphasized, their size 
is comparable to the experimental errors backing up our choice for the A^^. The nuclear 
modification for inclusive pion yield is defined as 



^dAu 



1 d^N^^^'/dprdy n^in^bias ^d'crj^'^/dpTdy 

d^al"^ /dprdy 



(iVeou) d^mVdpTdy 



where (iVcoii) denotes the number of binary nuclcon-nuclcon collisions and pt,^ the pion's 
transverse momentum and rapidity. A comparison with the PHENIX and STAR data is 



*See Appendix A of Ref. [7] for the detailed explanation. 
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Figure 2: The calculated i?^ and R^y compared with the NMC [UlIIl] and E772 [17] data. 
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Figure 3: Left: Comparison of the nuclear modifications Rq^ at = 1.69 GeV^ from 
HKN07 [5], nDS [6] and this work, EPS09 [7]. Right: The computed i?dAu for inclusive pion 
production compared with the PHENIX [18] and STAR [19] data multiplied by Jn = 1.03 
Jn = 0.90, respectively. 



shown in Fig. [3| Evidently, the shape of the spectrum — which in our calculation is a 
reflection of the similar shape in Rq — gets well reproduced by EPS09. Let us mention that 
the shape is practically independent of the fragmentation functions used in the calculation 
— modern sets like [T^[T31[I1] all give equal results. 

Figure [3] also presents a comparison of EPS09 gluon modifications Rq° with the earlier 
NLO analyses. The significant scatter of the curves highlight the difficulty of pinning down 
the nuclear modifications from the DIS and DY data alone — especially the behaviour of 
HKN07 looks different. Consequently, also the predictions for pion i?dAu differ significantly 
as is easily seen in Fig. |3j This is actually good news as this type of data, especially with 
better statistics, may eventually discriminate between different proposed gluon modifica- 
tions. 

Attention should be paid to the experimentally observed scaling-violations and to the 
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fact that the DGLAP dynamics reproduces them well. Most cleanly such effects are visible 
e.g from the small-x structure function ratios versus Q^, of which Fig. [¥] shows an example. 



1.0 

~; 0.95 
H 0.85 



f [GeV'] 1 



□ NMC 





EPS09NLO — 



.1=0.0175 

I I llll| HH 

10 



1.1 

- 1.05 

- 1.0 

- 0.95 

- 0.9 

- 0.85 

- 0.8 
— 1.1 



Figure 4: The calculated scale evolution of the ratio Ff"^ jF^ compared with the NMC data [20] . 

In summary, the very good agreement with the experimental data 'x^ /N w 0.79 found 
— especially the correct description of the scale-breaking effects — we argue, is compelling 
evidence for the applicability of coUinear factorization in nuclear environment. In addition to 
the best fit, we release [3T] 30 nPDF error-sets for practical use, encoding the neighborhood 
of the minimum. Although not discussed here, we have also performed the leading-order 
counterpart of the NLO analysis as we want to provide the uncertainty tools also for this 
widely-used framework. Although the best-fit quality is very similar both in LO and NLO, 
the uncertainty bands become smaller when going to higher order. 
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